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Abstract 


This  report  presents  trade-off  studies  on  Time-Frequency  Distribution  (TFD) 
algorithms  and  a  methodology  for  fusing  them  to  achieve  better  target  characterization. 
It  is  shown  that  TFD  algorithmic  fusion  considerably  increases  the  detectability  of 
signals  while  suppressing  artifacts  and  noise.  The  report  reviews  a  sample  of 
representative  TFD  algorithms.  Their  performance  is  studied  from  a  qualitative  and 
quantitative  point  of  view.  For  simplicity,  we  considered  the  mean-squared  error  as  a 
measure  of  performance  in  the  quantitative  trade-off  studies.  The  TFD  algorithmic 
fusion  is  performed  using  a  self-adaptive  signal.  It  may  be  adjusted  to  work  for  a  wide 
range  of  signal-to-noise  ratios.  The  algorithm  uses  the  first  two  terms  of  the  Volterra 
series  expansion  and  we  treat  the  outputs  of  the  time-frequency  algorithms  as  the 
variables  of  a  Volterra  series  and  the  coefficients  of  the  series  are  estimated  through 
training  sets  with  a  least-squares  algorithm.  Simplistic  TFD  algorithmic  fusion  methods 
(e.g.,  weighted  averaging  or  weighted  multiplication)  are  special  cases  of  the  proposed 
fusion  technique.  We  demonstrate  the  effectiveness  of  TFD  algorithmic  fusion  method 
using  experimental  High  Range  Resolution  (HRR)  radar  data  and  simulated  data. 
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I 


Resume 


Le  present  rapport  decrit  des  etudes  de  compromis  portant  sur  les  algorithmes  de 
distribution  temps-frequence  (DTF)  et  une  methode  permettant  de  les  fiisionner  en  vue 
d’accroitre  la  performance.  On  montre  que  la  fusion  d’algorithmes  de  DTF  accroit 
considerablement  la  detectabilite  des  signaux  tout  en  supprimant  les  artefacts  et  le  bruit. 
Le  rapport  examine  un  echantillon  d’algorithmes  de  DTF  representatifs.  Leur 
performance  est  etudiee  des  points  de  vue  qualitatif  et  quantitatif.  Par  souci  de 
simplicite,  nous  avons  utilise  l’erreur  quadratique  moyenne  comme  mesure  de  la 
performance  dans  les  etudes  de  compromis  quantitatives.  On  effectue  la  fusion  des 
algorithmes  de  DTF  en  utilisant  un  signal  auto-adaptatif.  Un  ajustement  peut  etre 
effectue  pour  F  utilisation  avec  une  vaste  gamme  rapports  signal/bruit.  En  utilisant  les 
deux  premiers  termes  du  developpement  en  serie  de  Volterra,  nous  traitons  les  resultats 
des  algorithmes  temps-frequence  comme  les  variables  d’une  serie  de  Volterra,  et  les 
coefficients  de  la  serie  sont  evalues  au  moyen  d’ensembles  d’entrainement  par  un 
algorithme  des  moindres  carres.  Les  methodes  simples  de  fusion  des  algorithmes  de 
DTF  (p.  ex.  moyennage  pondere  ou  multiplication  ponderee)  sont  des  cas  particuliers 
de  la  technique  de  fusion  proposee.  Nous  demontrons  l’efficacite  de  la  methode  de 
fusion  des  algorithmes  de  DTF  en  nous  servant  de  donnees  experimentales  du  radar  a 
haute  resolution  en  distance  (F1RR). 
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Executive  summary 


One  of  the  most  serious  deficiencies  in  NATO’s  air  defence  capability  and  the  one  that 
affects  almost  every  aspect  of  air  command  and  control,  and  weapon  systems  are  the 
lack  of  a  rapid  and  reliable  means  of  identifying  air  targets.  Radar  offers  a  capability  to 
identify  targets  over  long  distances,  under  all  weather  conditions,  during  day  or  night, 
without  the  need  for  the  target’s  cooperation.  This  mode  of  identification  is  known  as 
Non-Cooperative  Target  Recognition  (NCTR).  The  problem  of  NCTR  has  always  been 
a  topic  of  interest  in  military  operations,  that  is,  to  improve  situational  awareness. 

The  radar  echo  provides  a  target  profile  that  serves  as  a  ‘signature’  for  identification 
purposes.  There  are  a  number  of  radar  techniques  that  can  be  applied  to  NCTR.  Among 
some  of  the  more  promising  methods  are  High  Range  Resolution  (HRR),  Jet  Engine 
Modulation  (JEM)  and  Inverse  Synthetic  Aperture  Radar  (ISAR).  An  operational 
NCTR  system  must  satisfy  two  requirements:  accurate  identification  and  real  time 
operation.  It  has  been  identified  that  HRR  would  be  a  promising  candidate  for  an 
operational  NCTR  system.  HRR  has  a  relatively  simple  structure  (e.g.,  1-dimensional 
imagery)  and  has  an  all  aspect  capability.  It  requires  a  modest  signal-to-noise  ratio,  and 
is  applicable  to  a  large  class  of  radar  systems. 

Our  previous  experimental  studies  show  that  HRR  radar  image  profiles  can  be  severely 
distorted  when  a  target  possesses  very  small-perturbed  random  motions.  Therefore  the 
ability  to  generate  focussed  images  from  HRR  is  of  paramount  importance  to  military 
and  intelligence  operations.  One  of  the  central  problems  in  HRR  radar  data  is  the 
analysis  of  a  time  series.  Traditionally,  HRR  radar  signals  have  been  analysed  in  either 
the  time  or  the  frequency  domain.  The  Fourier  transform  is  at  the  heart  of  a  wide  range 
of  techniques  that  are  used  in  HRR  radar  data  analysis  and  processing.  The  Fourier 
transform-based  techniques  are  effective  as  long  as  the  frequency  contents  of  the  signal 
do  not  change  with  time.  However,  the  change  of  frequency  content  with  time  is  one  of 
the  main  features  we  observe  in  HRR  radar  data.  Because  of  this  change  of  frequency 
content  with  time,  HRR  radar  signals  belong  to  the  class  of  non-stationary  signals. 

C  onsequently,  for  the  inteipretation  of  radar  data  in  terms  of  a  changing  frequency 
content,  we  need  a  representation  of  our  data  as  a  function  of  both  time  and  frequency. 
Time-Frequency  Distributions  (TFDs)  describe  signals  in  term  of  their  joint 
time-frequency  content.  These  distributions  are  useful  for  analyzing  signals  with  both 
time  and  frequency  variations.  Therefore,  for  signals  with  time-varying  frequency 
contents,  TFDs  offer  a  powerful  analysis  tool.  Analysis  of  the  time-varying  Doppler 
signature  in  the  joint  time-frequency  domain  can  provide  useful  information  for  target 
detection,  classification  and  recognition. 

In  this  report,  a  new  TFD  algorithmic  fusion  method  has  been  presented  and  evaluated 
on  experimental  HRR  radar  data  and  simulated  data.  It  is  shown  that  our  TFD 
algorithmic  fusion  method  provides  an  effective  method  of  achieving  improved 
resolution,  highly  concentrated  and  readable  representation  without  the  auto-term 
distortion  and  cross-term  artifacts.  This  method  is  suitable  for  HRR  and  ISAR  data 
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where  multiple  scatterers  are  present,  and  that  noise  and  artifact  reduction  are  essential 
for  target  identification  applications.  Analysis  of  the  time-varying  Doppler  signature  in 
the  joint  time-frequency  domain  can  provide  useful  information  for  target  detection, 
classification  and  recognition.  We  anticipate  that  this  new  approach  will  find  a  wide 
range  of  uses  and  will  emerge  as  a  powerful  tool  for  time- varying  spectral  analysis. 

This  study  is  a  part  of  the  work  sponsored  by  the  US  Office  of  Naval  Research  to 
investigate  the  distortion  of  HRR  and  ISAR  images  under  the  US  Navy’s  International 
Collaborative  Opportunity  Program  (NICOP).  Canada  is  a  participant  in  this  project  on 
“Time-frequency  processing  for  ISAR  imaging  and  Non-Cooperative  Target 
Identification".  The  work  performed  in  this  report  is  relevant  to  the  ISAR  imaging 
capability  of  the  surveillance  radar  systems  on-board  of  the  CF  C-140  Aurora  patrol 
aircraft.  It  should  also  be  noted  that  target  recognition  based  on  radar  imagery  will  play 
an  important  role  in  future  CF  initiatives  on  ISR  (Intelligence,  Surveillance  and 
Reconnaissance)  for  land,  air  and  maritime  applications. 
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Sommaire 


Une  des  plus  graves  lacunes  de  la  capacite  de  defense  aerienne  de  l’OTAN,  qui  a  une 
incidence  sur  presque  tous  les  aspects  du  commandement  et  du  controle  aeriens  ainsi 
que  des  systemes  d’armes,  est  l’absence  d’un  moyen  rapide  et  fiable  d’identification  des 
objectifs  aeriens.  Le  radar  permet  d’identifier  les  objectifs  sur  de  grandes  distances, 
dans  toutes  les  conditions  meteorologiques,  le  jour  ou  la  nuit,  sans  necessiter  la 
cooperation  de  l’objectif.  Ce  mode  d’identification  est  appele  la  reconnaissance  de 
cible  non  cooperative  (NCTR).  Le  probleme  de  la  NCTR  a  toujours  presente  un  interet 
dans  les  operations  militaires,  par  exemple  du  point  de  vue  de  1’ amelioration  de  la 
connaissance  de  la  situation. 

L’echo  radar  foumit  un  profil  d’objectif  qui  sert  de  «  signature  »  aux  fins  de 
1’ identification.  II  existe  un  certain  nombre  de  techniques  radar  qui  peuvent  etre 
appliquees  a  la  NCTR.  La  technique  a  haute  resolution  en  distance  (HRR),  la 
modulation  des  moteurs  a  reaction  (JEM)  et  la  technique  du  radar  a  synthese 
d’ouverture  inverse  (ISAR)  comptent  parmi  les  methodes  les  plus  prometteuses.  Un 
systeme  de  NCTR  operationnel  doit  satisfaire  a  deux  exigences  :  identification  precise 
et  fonctionnement  en  temps  reel.  On  a  determine  que  la  technique  HRR  serait  une 
candidate  prometteuse  comme  systeme  de  NCTR  operationnel.  Elle  comporte  une 
structure  relativement  simple  (p.  ex.  imagerie  unidimensionnelle)  et  elle  possede  une 
capacite  tous  secteurs.  II  lui  sufiit  d’un  faible  rapport  signal/bruit  et  elle  est  applicable  a 
une  large  classe  de  systemes  radar. 

Nos  etudes  experimentales  anterieures  montrent  que  les  profils  d’images  obtenus  avec 
le  radar  HRR  peuvent  etre  fortement  deformes  lorsqu’un  objectif  effectue  des 
mouvements  aleatoires  tres  faiblement  perturbes.  Par  consequent,  la  capacite  de 
generer  des  images  mises  au  point  a  partir  du  radar  HRR  est  d’une  importance 
primordiale  pour  les  operations  militaires  et  les  operations  de  renseignement.  L’ analyse 
d’une  serie  chronologique  constitue  Tun  des  principaux  problemes  eprouves  avec  les 
donnees  du  radar  HRR.  Par  le  passe,  les  signaux  radar  ont  ete  analyses  soit  dans  le 
domaine  temporel,  soit  dans  le  domaine  frequentiel.  La  transformation  de  Fourier  est 
au  cceur  d’une  vaste  gamme  de  techniques  qui  sont  utilisees  pour  Tanalyse  et  le 
traitement  des  donnees  du  radar  HRR.  Les  techniques  a  transformation  de  Fourier  sont 
efficaces  dans  la  mesure  ou  le  contenu  frequentiel  du  signal  ne  varie  pas  avec  le  temps. 
Cependant,  la  variation  du  contenu  frequentiel  avec  le  temps  est  une  des  principales 
caracteristiques  que  nous  observons  dans  les  donnees  du  radar  HRR.  En  raison  de  cette 
variation  du  contenu  frequentiel  avec  le  temps,  les  signaux  du  radar  HRR  appartiennent 
a  la  categorie  des  signaux  non  stationnaires.  Par  consequent,  pour  T  interpretation  des 
donnees  radar  du  point  de  vue  d’un  contenu  frequentiel  changeant,  nous  avons  besoin 
d’une  representation  de  nos  donnees  en  fonction  du  temps  et  de  la  frequence.  Les 
distributions  temps-frequence(DTF)  decrivent  les  signaux  du  point  de  vue  de  leur 
contenu  temporel- frequentiel.  Elies  sont  utiles  pour  Tanalyse  des  signaux  qui 
presentent  des  variations  temporelles  et  frequentielles.  Par  consequent,  les  DTF 
constituent  un  outil  d’ analyse  puissant  des  signaux  dont  le  contenu  frequentiel  varie 


DRDC  Ottawa  TR  2003-187 


V 


dans  le  temps. 


Dans  ce  rapport,  nous  avons  presente  une  nouvelle  methode  de  fusion  des  algorithmes 
de  DTF  et  nous  l’avons  evaluee  en  l’appliquant  aux  donnees  experimentales  du  radar 
HRR.  On  montre  que  la  methode  de  fusion  des  algorithmes  de  DTF  constitue  un  moyen 
efficace  pour  obtenir  une  representation  a  resolution  accrue,  hautement  concentree  et 
lisible  sans  la  distorsion  des  auto-termes  ni  les  artefacts  de  termes  croises.  Cette 
methode  est  applicable  aux  donnees  du  radar  FIRR  et  du  radar  ISAR  lorsque  des 
diffuseurs  multiples  sont  presents,  et  que  la  reduction  du  bruit  et  des  artefacts  est 
essentielle  pour  les  applications  d’ identification  d’objectif.  Nous  esperons  que  cette 
nouvelle  approche  aura  aussi  des  applications  dans  une  vaste  gamme  de  situations  et 
qu’elle  se  revelera  un  outil  puissant  pour  l’analyse  des  spectres  variant  dans  le  temps. 
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1.  Introduction 


Our  previous  experimental  studies  show  that  High  Range  Resolution  (HRR)  radar 
image  profiles  can  be  severely  distorted  when  a  target  possesses  very  small-perturbed 
random  motions  [1],  The  ability  to  generate  focussed  images  from  HRR  is  of 
paramount  importance  to  military  and  intelligence  operations.  One  of  the  central 
problems  in  HRR  radar  data  is  the  analysis  of  a  time  series.  Traditionally,  HRR  radar 
signals  have  been  analysed  in  either  the  time  or  the  frequency  domain.  The  Fourier 
transform  is  at  the  heart  of  a  wide  range  of  techniques  that  are  used  in  HRR  radar  data 
analysis  and  processing.  The  Fourier  transform-based  techniques  are  effective  as  long 
as  the  frequency  contents  of  the  signal  do  not  change  with  time  [2-3].  However,  the 
change  of  frequency  content  with  time  is  one  of  the  main  features  we  observe  in  HRR 
radar  data  when  stepped  frequency  wave  form  HRR  processing  is  used.  Because  of  this 
change  of  frequency  content  with  time,  HRR  radar  signals  belong  to  the  class  of 
non-stationary  signals.  Consequently,  for  the  interpretation  of  radar  data  in  terms  of  a 
changing  frequency  content,  we  need  a  representation  of  our  data  as  a  function  of  both 
time  and  frequency.  Time-Frequency  Distributions  (TFDs)  describe  signals  in  term  of 
their  joint  time-frequency  content.  These  distributions  are  useful  for  analyzing  signals 
with  both  time  and  frequency  variations.  Therefore,  for  signals  with  time-varying 
frequency  contents,  TFDs  offer  a  powerful  analysis  tool. 

Over  the  past  ten  years,  radar  researchers  have  also  investigated  TFDs  as  a  unique  tool 
for  radar-specific  signal  analysis  and  image  processing  applications  [4-10].  It  was 
found  that  TFDs  provide  additional  insight  into  the  analysis,  interpretation,  and 
processing  of  radar  signals  that  are  sometimes  superior  to  what  is  achievable  in  the 
traditional  time  or  frequency  domain  alone.  The  specific  applications  where  TFDs  have 
been  used  include  signal  analysis  and  feature  extraction,  motion  compensation  and 
image  formation,  signal  SNR  (signal-to-noise  ratio)  improvement,  imaging  of  moving 
targets,  and  detection  of  moving  targets  [4-10]. 

The  most  widely  known  time-frequency  analysis  techniques  belong  to  the  Cohen  class 
with  the  leading  transform  being  the  Wigner-Ville  distribution,  which  is  based  on  a 
bilinear  model  [5-7,1 1-12].  This  model  is  able  to  analyze  time-varying  signals  with 
relatively  high  resolution.  However,  being  a  bilinear  model  it  introduces  cross-term 
artifacts[5-7,l  1-12],  Hence,  filtering  techniques  have  been  proposed  to  reduce  the 
cross-term  artifacts.  Among  them,  the  most  well  known  filtering  method  uses  the 
Choi- Williams  filter  [13].  Higher  than  bilinear  order  models  have  also  been  proposed 
for  time-frequency  analysis,  e.g.,  Wigner  bispectrum  and  trispectrum.  Choi-Williams 
filters  have  also  been  incoiporated  into  these  models.  These  higher-order  models  work 
well  with  isolated  signals  (time-varying  or  stationary).  However,  it  will  be  shown  in 
this  report  that  it  is  very  difficult  to  reduce  the  cross-term  artifacts  in  the  case  of 
multiple  signals  (i.e.,  multiple  point  scatterers  from  a  target). 

Wavelet-based  transforms  have  been  also  proposed  for  time-frequency  analysis. 
Wavelet  transforms  break  the  signal  into  a  set  of  bases  with  a  shape  that  is  based  on 
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affine  transformations  (i.e.,  translations  and  dilations)  of  a  basic  wavelet  called  the 
mother  wavelet.  Here  the  aspect  ratios  in  the  time-frequency  plane  are  proportional  to 
the  distance  from  the  zero  frequency  and  the  ratio  of  the  bandwidth  to  the  center 
frequency  remains  constant  [4],  Among  the  most  simple  but  also  effective 
time-frequency  wavelets  is  the  Gabor  based  time-frequency  expansion,  where  the 
mother  wavelet  is  a  Gaussian  function.  The  Gabor  transform  originally  used  rectangles 
to  designate  each  of  the  time- frequency  elementary  signals.  Wavelet-based 
time-frequency  analysis  belongs  to  the  Weyl-Hiesenberg  generalized  class.  The 
introduction  of  chirplets  for  time-frequency  analysis  also  belongs  to  the  same  class  and 
has  been  used  to  better  describe  the  second-  (i.e.  Doppler  acceleration)  and 
higher-order  signal  nonlinearity  [5,14-15].  Another  transform  that  is  also  considered  as 
a  generalization  of  the  Gabor  transform  is  the  s-transform  [16-17].  All  these  transforms 
do  not  introduce  cross-term  artifacts  as  those  of  the  Cohen  class  time-frequency 
analysis  methods.  However,  they  often  result  in  inferior  resolution  when  compared  to 
the  Cohen  class  transforms. 

Other  techniques  in  time-frequency  analysis  include  the  Adaptive  Joint  Time-Frequency 
(AJTF)  transform  [8-9,18-19].  It  is  based  on  the  series  analysis  of  a  time-varying  signal 
into  a  high-order  polynomial  of  time  (i.e.,  Doppler,  Doppler  acceleration,  etc.).  The 
estimation  of  parameters  assumes  the  presence  of  very  strong  (i.e.,  high  signal-to-noise 
ratio)  point  scatterers.  For  example,  for  a  second-order  polynomial,  at  least  two  very 
strong  scatterers  are  required  from  the  same  moving  target.  These  polynomial 
parameters  may  be  used  for  motion  compensation  as  well  as  for  target  recognition. 

In  recent  years  the  Fractional  Fourier  Transform  (FRFT)  has  been  introduced.  The 
FRFT  is  able  to  find  linear  changes  in  the  frequency  overtime  [8].  However,  for 
higher-order  nonlinearity,  extensions  of  this  transform  to  time-frequency  domain  must 
be  incorporated  [14], 

Enhancements  and  filtering  techniques  on  the  results  of  the  above  mentioned 
transforms  have  been  proposed  by  many  researchers.  The  most  promising  methods  are 
the  Choi- Williams  filters  which  are  applied  in  the  Wigner-Ville  class  of  time-frequency 
transforms  and  the  reassignment  methods,  which  are  applied  to  any  time- frequency 
distribution  [20-21]. 

In  this  report  we  present  a  new  method  for  increasing  the  detectability  of  time-varying 
signals,  for  high  resolution  improvements  on  the  results  and  for  reduction  of  the  noise 
and  cross-term  artifacts.  This  new  concept  of  time-frequency  analysis  is  based  on  TFD 
algorithmic  fusion.  TFD  algorithmic  fusion  is  required  to  mathematically  analyze  the 
TFD  algorithms.  This  analysis  may  be  performed  a  priori,  or,  by  using  self-adaptive 
data  analysis  techniques  on  the  results  (e.g.,  by  applying  weighting  techniques,  where 
the  weights  are  found  optimally  from  training  data  sets  using  least-squares  techniques 
or  neural  networks). 
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2.  Time-Frequency  Analysis  Algorithms 


In  this  section  we  present  a  sample  of  algorithms  used  in  our  performance  evaluation 
and  the  development  of  an  algorithmic  fusion  method  to  enhance  the  results.  Detailed 
description,  properties  and  derivation  of  each  transform  can  be  found  from  the 
references,  particularly  from  [4-13]. 

2.1  Adaptive  Energy  Distribution  (AED) 

Adaptive  Energy  Distribution  (  AED)  is  a  window-matched  spectrogram,  where  the 
window  matching  is  performed  with  the  use  of  generalized  instantaneous  parameters 
(GIP’s).  The  GIP’s  that  are  used  for  Gaussian  window  matching  are  the  instantaneous 
time  variance  (ITV),  the  instantaneous  frequency  variance  (IFV),  and  the  instantaneous 
time-frequency  (ITF)  covariance  [22].  These  variance  measures  are  calculated  from 
instantaneous  moments  and  expressed  in  the  following  succinct  forms,  i.e., 


(1)  ITVzh(t,  f) 


(2)  IFVzh(t,f) 


1 

87 r2 


5ft 


’ d2Sx(t,f ) 

1 

(dSx(t,f) 

1  VI 

dp 

sptj) 

V  df 

~ d2Sx(t,  f) 

1 

f  dSx(t,  f) 

1  yl 

dt 2  Sx(t,f)  \  dt  Sx(t,f)J 


d2sx(tj )  1  _  dsx(tj)  dSx(t,  f)  f  1  Y 

dfdt  Sx(t,f)  df  dt  \Sx(t,f)J 
(3) 

where  Sx(t,  /)  specifies  the  spectrogram  of  the  signal  x(t)  formed  using  the  window 
h(t).  5ft  indicates  the  real  part  of  the  expression  and  /  =  u/ {2i r),  where  /  is  the 
frequency  and  c o  is  the  angular  frequency.  The  subscript  zh  in  equations  (l)-(3)  denotes 
that  the  signal  x(t)  has  been  multiplied  by  the  window  hit).  The  parameters  for  an 
optimal  window  are  found  by  matching  the  GIP’s  of  the  signal-window  spectrogram 
to  those  of  the  window-window  spectrogram.  The  optimal  window  parameters  are 
updated  iteratively  until  the  desired  level  of  matching  is  achieved.  At  the  nth  iteration 
they  are  calculated  with 


ITFzh(t,f)  =  5ft 


(4) 

and 


an+ 1  (t,  /) 


(Svr)3  (, TTFzhn(t,f))2IFVzhn(t,f ) 
(87t)2  (ITFzhn(t,  f))2  +  1 
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(5) 


(87r)2  ITFzhn  (t,  f)  IFVzhn  (t,  f) 


Pn+l(t,f) 


(8*)2(ITFehn(t,f))*  +  l 


or  if  ITVzhn  (t,  f )  >  IFVzhn  (t,  f)  then 


(6) 

and 


an+ 1  (t,  /) 


1 

8t r-ITVzhn(t,f) 


(7) 


Pn+l(tJ) 


1 

(87t)2  ITVzhn  (t,  f)  ITFzhn  (t,  f) 


where  a  and  / 3  are  the  parameters  of  the  Gaussian  window  h(t)  =  exp 


n(a-j/3)t2 


2.2  Adaptive  Gaussian  Chirplet  Decomposition  (AGCD) 

Adaptive  joint  time-frequency  analysis  is  a  good  tool  to  analyze  the  time-varying 
Doppler  features  of  a  target.  The  chirp  function  is  one  of  the  most  fundamental 
functions  in  nature.  Many  natural  events,  for  example,  signals  in  seismology  and  radar 
systems,  can  be  modelled  as  a  superposition  of  short-lived  chiip  functions.  Hence,  the 
chiip-based  signal  representation,  such  as  the  Gaussian  chirplet  decomposition,  has 
been  an  active  research  area  in  the  field  of  signal  processing  [4,15].  The  Gaussian 
chirplet  is  defined  as 


(8)  hk(t)  =  exp  |-^  (t  -  tkf  +j  ^  +  y(t-  4)^  (t  -  4) | 

where  ak  >  0  and  4>  ^k^Pk£^-  hk{t)  has  a  very  short  and  smooth  envelope,  i.e.,  a 
Gaussian  envelope.  (4,  wk)  indicates  the  time  and  frequency  center  of  the  linear  chirp 
function.  The  variance  ak  controls  the  width  of  the  chirp  function.  The  parameter  j3k 
determines  the  rate  of  change  of  frequency.  It  shows  that  not  only  can  we  adjust  the 
variance  and  time-frequency  center,  but  also  we  can  regulate  the  orientation  of  hk{t)  in 
the  time-frequency  domain  by  varying  parameter  (3k ■  F°r  the  Gaussian  chirplet 
decomposition,  the  corresponding  adaptive  spectrogram  is 


(9) 


ASM*  2]T 
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\Bk\2  exp 


~ak  {t  -  4)2 


—  (w  -  uk  ~  (3kt)2 

Oik 


which  is  non-negative.  AGCD  has  an  excellent  resolution  and  has  no  cross-term 
artifacts. 
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2.3  Adaptive  Gabor  Expansion  (AGE) 


Adaptive  Gabor  Expansion  (AGE)  is  similar  to  the  adaptive  Gaussian  chirplet 
decomposition  [4,15].  The  only  difference  is  AGE  uses  a  Gaussian  function  instead  of  a 
Gaussian  chirplet  function,  i.e., 


(10)  hk{t) 


{~~Y  (*  “  4)2  +  juk  (t  -  4)|  tk,ujkeR,akeR+ 


2.4  Linear  Time-Frequency  Transforms 

All  linear  time-frequency  transforms  satisfy  the  supeiposition  or  linearity  principle 
which  states  that  if  x(t)  is  a  linear  combination  of  some  signal  components,  then  the 
time-frequency  transforms  of  x(t)  have  the  same  linear  combination  of  the 
time-frequency  transforms  of  each  of  the  signal  components,  i.e., 


(11)  x(t)  =  C!  Xi(t)  +  C2  x2(t)  -*■  Tx(t,f)  =  Cl  TX1  (t,  f)  +c2TX2(t,f ) 

Linearity  is  a  desirable  property  in  any  application  involving  multi  component  signals. 
One  of  the  basic  linear  time-frequency  transforms  is  the  short-time  Fourier  transform 
(STFT)  [5-7,10].  The  STFT  is  given  as 


(12) 


STFTx(t,co) 


x(r)  w(t  —  t )  exp(— jcor)  dr 


where  w(t)  is  a  rectangular  window  function,  t  is  time  and  r  is  running  time.  When  a 
Gaussian  window  is  used,  the  result  is  the  Gabor  transform  and  is  given  as 


(13)  GTFTx(t,uj)  =  x(t)  --L=  exp  |  |  exp(-jo;r)  dr 

where  a  is  the  standard  deviation.  Other  windowing  functions  that  are  used  in  this 
report  are  [5-7]: 


(14) 


Triangular  window  :  w(t) 
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1*1  < 
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(15) 


Hanning  window  :  w(t) 
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(16) 


Hamming  window  : 


w(t)  =  0.54  +  0.46  cos 


2irt 


for 


1*1  < 


T 

2’ 


(17) 


Blackman  window  : 


2yrt  2tt) t 

w(t )  =  0.42—0.5  cos - 1-0.08  cos - 

r  r 


for 


In 


(18) 


Kaiser  window  :  w(t)  = 


7V7!  -  (2t/r)2 


4(7) 


for  |f|  < 


2’ 


where  1 0  is  the  modified  Bessel  function  of  the  first  kind  and  of  order  zero,  and  7  is  a 
parameter.  We  also  introduce  and  test  an  inverse  relation  window  function.  Its 
amplitude  response  relatively  has  a  narrow  main- lobe  and  has  almost  no  side-lobes. 
The  inverse  window  function  is  given  as 


(19) 


Inverse  window  function  :  w(t) 


where  the  exponential  p  is  a  positive  integer. 


2.5  Bilinear  Time-Frequency  Distributions 

In  contrast  to  the  linear  time-frequency  transforms  such  as  the  short-time  Fourier 
transform,  the  Wigner-Ville  distribution  (WVD)  is  said  to  be  bilinear.  For  a  signal  s(t), 
its  Wigner-Ville  distribution  is 


(20)  Wx(t,  js)  -  I  x(t  +  ^)  x*(t  -  ^)  exp(-jwr)  dr 

The  Wigner-Ville  distribution  not  only  possesses  many  useful  properties,  but  also  has 
better  resolution  than  the  STFT  spectrogram.  Although  the  Wigner-Ville  distribution 
has  existed  for  a  long  time,  its  applications  are  very  limited.  One  main  deficiency  of  the 
Wigner-Ville  distribution  is  the  so-called  cross-term  artifacts  [5-7].  Suppose  we  express 
a  signal  as  the  sum  of  two  signal  components, 


(21) 


x(t)  =  Xl  (t)  +  X2{t) 
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Substituting  this  into  equation  (20),  we  have 


(22) 


Wx  (^5  k'O  —  ^^X\X\  (^5^)  “1“  ^^X2X2  (^5  ”1”  ^#1#2  (^5  ^"0  "I” 


where 


(23)  n'XlX2(/.a,')  -  ^  J  xl(t  -  ^t)  X2(t  +  ^t)  e  JUJT  dr 

This  is  called  the  cross  Wigner-Ville  distribution.  In  terms  of  the  spectrum  it  is 

(24)  WxlX2(t,u)  =  ^  J  x\{u  +  ^6)x2(u  -  ^9)e~jted6 

The  cross  Wigner-Ville  distribution  is  complex.  However,  WXlX2  =  W*2X1,  and 
therefore  WX1X2  (t,  u>)  +  WX2X1  (t,  tu)  is  real.  Hence 


(25)  W  (t,  to)  =  WXlXl  ( t ,  uj)  +  WX2X2  ( t ,  io)  +  2  3?{WXlX2  (f ,  w)} 

We  see  that  the  Wigner-Ville  distribution  of  the  sum  of  the  two  signals  is  not  the  sum  of 
the  Wigner-Ville  distribution  of  each  signal  but  has  the  additional  term 
2  $i{ WXlX2  (t.  tu)}.  The  term  is  often  called  the  interference  term  or  the  cross-term  and 
it  is  often  said  to  give  rise  to  artifacts.  Because  the  cross-term  usually  oscillates  and  its 
magnitude  is  twice  as  large  as  that  of  the  auto-terms,  it  often  obscures  the  useful 
time-dependent  spectrum  patterns. 

The  inverse  Fourier  transform  of  the  Wigner-Ville  distribution  is  called  the  ambiguity 
function  (AF).  The  Fourier  transform  maps  the  Wigner-Ville  distribution  auto-terms  to 
a  region  centered  on  the  region  of  the  AF  plane,  whereas  it  maps  the  oscillatory 
Wigner-Ville  distribution  cross-terms  away  from  the  origin  [5-7,23-24]. 

The  fact  that  the  auto-  and  cross-terms  are  spatially  separated  in  the  AF  domain  means 
that  if  we  apply  a  filter  function  to  the  AF,  we  can  suppress  some  of  the  cross-terms. 
This  filtering  operation  defines  a  new  time- frequency  distribution 


(26)  TFD  =  Fourier  transform! AF  •  Kernel} 

with  properties  different  from  the  Wigner-Ville  distribution.  The  filter  function  is  called 
the  ‘kernel’  of  the  TFD.  Since  there  are  many  possible  2-dimensional  kernel  functions, 
there  exist  many  different  TFDs  for  the  same  signal.  The  class  of  all  TFDs  obtained  in 
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this  fashion  is  called  Cohen’s  class.  A  more  detailed  description  of  the  ambiguity 
function  and  Cohen’s  class  are  given  [5-7,23-24], 

Cohen’s  general  class  of  bilinear  TFDs  is  defined  as 


(27) 


C(t,u) 


AF(6,t)  <!>(#,  t)  exp \j  (0t  —  ujt)]  d0  dr 


where  <J>(0,  r)  is  a  kernel  function  and  AF(0.  r)  is  the  ambiguity  function  of  the  signal 
defined  by 


(28)  AF(0,t)  =  J  x(t  +  7^)x*{t  -  ^)  exp(-j0T)dt 

where  6  represents  the  frequency-offset  and  r  represents  the  time-offset.  Cross-terms  in 
the  WVD  can  be  effectively  filtered  in  the  ambiguity  domain  by  designing  kernels  that 
filter  the  signal  part  of  the  ambiguity  function.  The  kernel  functions  that  are  used  in  this 
report  are  [5-7]: 

Bom-Jordan: 


(29) 

Choi- Williams: 


$(0,t)  =  2 


sin(#r/2) 

0T 


(30) 


$(#,  r)  =  exp 


Butterworth: 


(31) 


$(<9,t) 
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Generalized  exponential: 


(32) 


<h(0,  r)  =  exp 


2  TV "1  ~1 


Cone-shaped, 


8 


DRDC  Ottawa  TR  2003-187 


(33) 


4>(0,t)  =  2  exp(— ar2) 


sin(0r/2) 

0 


Tilted  Gaussian: 


(34) 


<F(0,  r)  =  exp 


Tilted  Butterworth 


(35) 


$(0,t) 


t9 

tq00 


-l 


where  r o  and  do  are  two  scaling  parameters  and  r  is  the  rotation  parameter.  The 
kernel’s  contours  are  un-tilted  when  r  =  0  and  tilted  ellipses  when  0  in  the  AF 
plane.  A  tilted  parallel  strip  is  also  possible  by  setting  r  ±  1. 


2.6  Adaptive  Optimal-Kernel  (AOK) 

Since  real  signals  have  different  shapes  in  the  ambiguity  domain,  no  single  kernel  can 
give  adequate  performance  for  a  large  class  of  signals.  Flence,  there  has  been  increasing 
interest  in  signal-dependent  or  adaptive  TFDs,  in  which  the  kernel  function  varies  with 
the  signal.  Adaptation  of  the  kernel  over  time  is  beneficial  because  it  permits  the  kernel 
to  match  the  local  signal  characteristics  [25]. 

The  signal-dependent  TFR  is  based  on  kernels  with  Gaussian  radial  cross-sections 

(  02  -f-  T2  \ 

(30)  mT>=exp{-2^)) 

The  function  cr(ip)  controls  the  spread  of  the  Gaussian  function  at  radial  angle  ip. 
Clearly  if  a  is  smooth,  then  4>  is  smooth.  The  angle  ip  is  measured  between  the  radial 
line  through  the  point  ( 6 ,  r)  and  the  0-axis 


(37) 


ip  =  arctan 


6 


It  is  natural  to  express  radially  Gaussian  kernels  in  polar  coordinates;  using 
r  =  y^T  t2  as  the  radius  variable, 
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(38) 


<b(r,  ip)  =  exp 


2  a2  (ip) 


A  better  TFD  results  when  the  kernel  is  well  matched  to  the  components  of  a  given 
signal.  The  radially  Gaussian  kernel  is  adapted  to  a  particular  signal  by  solving  the 
following  optimization  problem: 


/*27T  /‘OO 

(39)  max  /  /  \AF(r,  ip)  4>(r,  ip)\2  r  dr  dip 

$  Jo  Jo 

Subject  to 

(40)  <f>(r,  ip)  =  exp  >  J  0-2  WO  dW  <  «  >  0 

Here,  AF(r,  ip)  is  the  AF  of  the  signal  in  polar  coordinates.  The  problem  requires  the 
area  of  the  kernel  to  be  constant  with  the  projection  of  the  kernel  on  the  ambiguity 
function  of  the  signal  to  be  maximal.  Since  this  optimization  problem  takes  a  long  time 
to  solve,  a  simplified  adaptive  kernel  is  implemented.  The  simplified  adaptive  kernel  is 
just  proportional  to  the  sum  of  the  ambiguity  function  in  each  radial  direction,  i.e., 


POO 

(41)  a(ip)  =  c  /  \AF(r,  ip)\2  dr,  c>  0 

Jo 

The  parameter  c  is  similar  to  a  in  equation  (40)  and  determines  the  extent  of  the  kernel 
in  a  radial  direction.  This  way  the  kernel  lets  through  most  of  the  signal’s  auto-terms 
but  not  the  cross-terms  and  calculates  the  optimization  problem  more  quickly.  An 
adaptive  kernel,  which  is  radially  inverse,  has  been  also  proposed  for  the  TFD,  i.e., 


(42) 


4>(r,  ip) 


1  + 


2 


-1 


An  adaptive  kernel,  which  is  radially  constant,  is  also  proposed  for  the  TFD  and  is 
given  as 


(43) 


<F(r,V’) 


1  if  expf-2^))^fc 

0  if  exp(-2^))<fc 


where  k  is  a  constant 
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3.  Fusion  of  algorithms  -  A  new  approach 


Fusion  of  TFD  algorithms  requires  mathematical  analysis  of  the  performance  of  each 
algorithm  under  a  wide  range  of  single  or  multiple  signals  of  interest  at  all 
signal-to-noise  ratios  and  with  noise  models  that  are  representative  to  the  operating 
environment  of  applicability.  The  objective  of  fusing  TFD  algorithms  is  to  develop  a 
system  that  has  better  performance  than  any  subset  of  algorithms  or  their  combinations. 

It  is  extremely  difficult  to  develop  a  single  TFD  algorithm  that  produces  satisfactory 
resolution  and  low  noise  levels  with  no  cross-term  artifacts  for  all  signals  and 
signal-to-noise  ratios.  However,  different  TFD  may  filter  various  components  of  signals 
with  different  degree  of  success,  based  on  the  quality  and  number  of  signals  and  their 
signal-to-noise  ratios.  With  many  TFD  algorithms  available  a  combination  of  them  may 
result  in  both  high  resolution  and  low  levels  of  noise  and  cross-term  artifacts.  TFD 
algorithmic  fusion  is  designed  to  amplify  results  from  algorithms  that  produce  good 
TFDs  and  to  suppress  algorithms  that  do  not  work  well  for  a  particular  environment. 

Mathematical  analysis  of  the  TFD  algorithms,  for  all  kinds  of  environments,  is  not  an 
easy  task.  In  this  section  we  present  a  methodology  that  weights  the  outputs  of  selected 
TFDs  and  their  first-order  multiplication  (i.e.  first-  and  second-order  terms  of  a  Volterra 
series  expansion)  to  result  in  an  enhanced  TFD.  The  enhanced  TFD  produces  high 
resolution,  low  cross-term  artifacts  and  low  noise.  Herein,  the  presented  fusion  process 
produces  a  new  TFD,  which  is  a  linear  combination  of  the  TFDs  from  available 
algorithms  and  their  products. 

The  fusion  process  finds  the  coefficients  pi  for  each  of  the  TFDs  A/  =  { }  /  of  the 
available  algorithms  such  that  the  linear  combination  is  a  close  as  possible  to  the  ideal, 
noise  free,  TFD: 


(44) 


P\A\  +  P2A2  +  ...  +  pnAn  —  D 


Equivalently  equation  (44)  can  be  written  as 


(45) 


^  '  Pl{aij)l  —  dij 

1=1... n 


where  D  =  {(%}  is  the  ideal  distribution.  In  this  description,  n  refers  to  the  number  of 
available  time- frequency  algorithms  and  their  TFDs,  and  N  is  the  length  of  a  signal.  To 
solve  for  p  each  Ai  matrix  is  written  as  a  column  and  the  equation  becomes 
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(46) 
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which  has  the  least-square  solution 

(47)  p  =  [BtB]~1Btc 

where  B  is  the  matrix  that  contains  the  TFDs  of  each  algorithm  as  columns  and  c  is  the 
column  vector  that  contains  all  the  dij.  Once  all  the  coefficients  p/  are  known,  a  fused 
TFD  F  can  be  constructed  according  to 


(48)  Pl-A-l  +  P2-A-2  +  •••  +  Pn-A-n  —  F 

To  include  multiplications  between  the  TFDs  of  different  algorithms,  equation  (45)  can 
be  changed  to 


(49)  ^  piAi  +  PijAiAj  —  D 

1=1... n  i=l...rrJ=i-\-l...n 

where  AiAj  =  A  are  element  by  element  products  between  each  pair  of  TFDs  and  pk 
are  their  coefficients  in  the  fused  TFD.  With  this  modification,  the  matrix  B  will 
include  the  columns  of  Ai  as  well  as  /l/,.  and  p  will  be  a  column  vector  with  both  pl  and 
pk  in  it.  Solution  of  equation  (47)  will  give  the  coefficients  that  are  required  to  make 
the  fused  TFD.  When  multiplicative  factors  are  included  in  the  fusion  process,  the 
signal  parts  are  amplified,  and  at  the  same  time,  blurring  and  artifacts  are  suppressed. 

In  order  to  make  the  fusion  process  to  suppress  noise  (in  addition  to  artifacts)  in  the 
TFDs  of  the  given  algorithms,  the  coefficients  p  have  to  include  the  noise  information. 
When  TFDs  with  different  noise  levels  (SNR)  are  present,  then  a  separate  matrix 
gSNRm  js  forme(j  for  each  SNR,  i.e., 
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Then  equation  (46)  becomes 


(51) 


'  B3NR1  ' 

bsnr2 

c 

c 

P  = 

BSNRm 

c 

Solution  of  equation  (5 1)  can  then  be  used  to  make  a  TFD  in  which  noise  and  artifacts 
are  suppressed.  Equation  (51)  may  be  extended  to  include  multiple  signals  and 
furthermore  make  the  proposed  approach  applicable  to  a  wider  set  of  signals.  It  should 
be  pointed  out  that  equation  (51)  is  used  for  the  training  process,  where  only  the  vector 
p  is  unknown.  The  left  part  of  the  equation  (51)  is  used  for  the  fusion  process,  when  the 
signal  is  unknown. 


Although  the  fusion  process  requires  an  ideal  TFD  in  order  to  find  the  coefficients  of 
each  available  distribution,  it  is  required  only  once.  After  the  coefficients  are  estimated, 
they  are  used  to  construct  fused  TFDs  of  signals  from  the  same  source,  but  with 
different  levels  of  noise  (variable  signal-to-noise  ratios)  and  a  set  of  representative 
noise  models  and  signals. 
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Figure  1:  (a)  Time  history  of  bat  echoiocation  signal,  (b)  Fourier  spectrum  of  bat  echolocation  signal. 

4.  Results  and  discussion 

4.1  Performance  evaluation  on  time-frequency  algorithms 

In  this  section  we  evaluate  and  compare  the  qualitative  and  quantitative  performances 
of  different  TFDs  with  an  ideal  distribution.  A  TFD  maps  a  one-dimensional  signal  to  a 
two-dimensional  time- frequency  image  that  displays  how  the  frequency  content  of  the 
signal  changes  over  time.  Experimental  bat  echolocation  signals  provide  an  excellent 
motivation  for  time-frequency  based  signal  processing.  The  bats  use  echolocation 
signals  for  navigation  and  hunting  in  a  dark  environment.  Bat  echolocation  signals  are 
high  frequency  and  short  duration  signals.  These  signals  are  made  up  of  chirps  having 
hyperbolic-like  instantaneous  frequency.  For  this  signal,  the  frequencies  tend  to 
increase  with  increasing  time.  The  experimental  echolocation  signal  data  was  taken 
from  [26]  and  is  used  for  the  comparison  and  evaluation  of  different  TFDs.  There  are 
400  data  samples  with  sampling  period  of  7  fis.  Therefore  the  data  length  is  2.8  ms.  In 
this  section  we  evaluate  the  ability  of  each  TFD  to  detect  the  hyperbolic  nature  of  the 
chirp  with  minimum  artifacts. 

Figure  1  shows  that  neither  the  time  signal  nor  its  Fourier  spectrum  reveal  the  true 
structure  of  the  signal.  In  contrast,  a  time-frequency  image  of  the  signal  clearly  exposes 
its  non-stationary  character.  While  each  signal  has  a  unique  Fourier  spectrum,  a 
time-frequency  analysis  of  that  signal  is  not  unique.  In  other  words,  many  different 
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Figure  2:  An  ideal  time-frequency  distribution.  The  horizontal  axis  is  time  and  the  vertical  axis  is  frequency. 


TFDs  can  describe  the  same  data.  Since  for  any  given  signal  some  TFDs  are  ‘better 
than  others’,  TFD  design  has  become  an  important  research  area. 

An  ideal  distribution  or  the  theoretical  instantaneous  frequency,  i.e.,  noise-free 
Wigner-Ville  distribution,  of  a  bat  echolocation  signal  is  manually  plotted  in  Figure  2  as 
a  reference  for  comparison.  The  qualitative  performance  of  different  TFDs  was 
examined  based  on  the  amount  of  cross-term  artifacts,  noise  levels  and  time-frequency 
resolution  (Figures  3-6).  The  time- frequency  resolution  of  different  TFDs  was  visually 
compared  to  the  ideal  time-frequency  distribution.  Although  we  studied  and  compared 
all  the  TFDs  described  in  section  2,  we  only  show  certain  TFDs  in  Figures  3-6. 
However,  the  quantitative  performances  of  all  the  TFDs  are  tabulated  in  Tables  1-2.  The 
quantitative  performance  of  each  TFD  was  studied  using  the  mean-squared  error 
estimation.  The  mean-squared  error  was  estimated  between  TFD  produced  by  each 
algorithm  and  an  ideal  distribution. 


Figure  3  shows  the  adaptive  energy  distribution  and  adaptive  spectrograms  computed 
by  Gabor  expansion  and  Gaussian  chirplet  decomposition  of  a  bat  signal.  In  this 
example  the  time-frequency  resolution  of  the  adaptive  Gaussian  chirplet  spectrogram  is 
obviously  superior  to  other  two  schemes.  The  adaptive  Gaussian  chirplet  spectrogram 
not  only  shows  excellent  time-frequency  resolution,  but  also  does  not  suffer  from 
cross-term  artifacts.  Compared  with  the  adaptive  energy  distribution,  the  performance 
of  the  adaptive  Gabor  expansion  based  spectrogram  is  rather  poor.  In  this  case,  adaptive 
Gabor  expansion  based  spectrogram  does  not  offer  good  time-frequency  resolution  due 
to  the  limitation  of  the  number  of  degrees  of  freedom  used  in  the  elementary  functions. 
Compared  with  the  adaptive  Gabor  transform,  the  adaptive  chirplet  matches  the  bat 
sound  much  better  in  terms  of  fewer  terms  required.  In  this  case,  the  adaptive  Gabor 
transform  needs  1 0  times  more  terms  than  that  needed  by  its  adaptive  Gaussian  chirplet 
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Figure  3:  Time-frequency  distribution  of  (a)  AED,  (b)  AGCD  and  (c)  AGE.  The  horizontal  axis  is  time  and  the  vertical 
axis  is  frequency. 


counterpart.  This  is  because  the  Gaussian  chirplet  fits  the  bat  sound  better  than  it  does 
the  regular  Gaussian  function.  The  regular  Gaussian  function  is  a  special  case  of  the 
Gaussian  chirplet.  The  adaptive  Gaussian  chirplet  decomposition  algorithm  is  an 
excellent  candidate  for  signals  that  are  composed  of  chirps  and  chirplets,  but  for  signals 
that  have  frequency  modulation,  the  adaptive  Gaussian  chirplet  decomposition  is  not 
the  desired  algorithm  since  it  tries  to  approximate  the  signal  by  chirplets  (line  segments 
in  the  TF  plane,  see  Figure  3b). 

Figure  4  shows  the  linear  time-frequency  transforms  with  different  windows.  In  this 
case  short-time  Fourier  transform  is  used  with  various  window  functions.  The  window 
functions  used  in  this  figure  are  Gaussian,  Hanning  and  Kaiser.  Unlike  the  adaptive 
Gaussian  chirplet  decomposition,  these  transforms  represent  all  the  information 
contained  in  the  signal.  These  linear  time-frequency  transforms  often  work  well,  but  the 
performance  always  depends  on  the  choice  of  the  window  used  in  the  representation.  In 
general,  the  appropriate  window  depends  on  the  data  and  can  differ  for  different 
components  in  the  same  signal.  Furthermore,  selection  of  the  appropriate  window 
requires  some  knowledge  of  the  signal  components  of  interest.  In  many  applications, 
such  information  is  not  available,  and  it  is  desirable  to  avoid  presupposing  the  form  of 
the  data.  This  is  one  of  the  main  drawbacks  of  linear  time-frequency  transforms. 

Figure  5  shows  the  bilinear  TFDs  with  different  kernels.  The  kernel  functions  used  in 
this  figure  are  Choi- Williams,  tilted  Butterworth  and  tilted  Gaussian.  These  bilinear 
TFDs  generally  show  better  time-frequency  resolution  than  the  linear  transforms.  But 
these  bilinear  distributions  require  longer  computational  time  and  have  larger  amounts 
of  cross-term  artifacts.  Since  real  signals  have  different  ambiguity  functions,  no  single 
kernel  can  perfectly  filter  the  cross-term  artifacts.  When  a  low-pass  kernel  is  employed, 
there  is  a  trade-off  between  cross-term  suppression  and  auto-term  concentration. 
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Figure  4:  Linear  time-frequency  distributions  with  different  windows:  (a)  Gaussian,  (b)  Hanning  and  (c)  Kaiser.  The 
horizontal  axis  is  time  and  the  vertical  axis  is  frequency. 


Figure  5:  Bilinear  time-frequency  distributions  with  different  kernels:  (a)  Choi-Williams,  (b)  tilted  butterworth  and  (c) 
tilted  Gaussian.  The  horizontal  axis  is  time  and  the  vertical  axis  is  frequency. 
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Figure  6:  Bilinear  time-frequency  distributions  with  adaptive  kernels:  (a)  radially  constant,  (b)  radially  Gaussian  and 
(c)  radially  inverse.  The  horizontal  axis  is  time  and  the  vertical  axis  is  frequency. 


Generally,  as  the  pass-band  region  of  the  kernel  reduces,  the  amount  of  cross-term 
artifacts  suppression  increases,  but  at  the  expense  of  auto-term  concentration. 
Adjustable  kernels  such  as  tilted  Gaussian  and  tilted  Butterworth  can  be  manually  fitted 
to  a  signal  for  better  results.  However,  for  some  signals,  the  selection  of  manual 
parameters  to  produce  good  results  is  very  difficult.  Figure  5  shows  that  the  tilted 
Gaussian  kernel  function  produces  the  best  results  in  terms  of  low  cross-term  artifacts 
and  relatively  higher  time-frequency  resolution  than  Figures  3-4. 

In  order  to  compare  the  performance  of  the  signal-dependent  TFD  with  the  fixed-kernel 
bilinear  TFDs  shown  in  Figure  5,  the  adaptive  optimal  kernel  TFDs  are  shown  in  Figure 
6.  The  kernel  functions  used  in  this  figure  are  radially  constant,  radially  Gaussian  and 
radially  inverse.  Unlike  the  kernels  used  in  bilinear  TFD  that  emphasize  preserving  the 
properties  of  the  Wigner-ville  distribution  over  matching  the  shape  of  auto-terms,  the 
signal-dependent  kernels  aim  to  optimally  pass  the  auto-terms  while  suppressing 
cross-term  artifacts.  Since  a  fixed-kernel  acts  on  the  ambiguity  domain  as  a  filter,  it  is 
limited  in  its  ability  to  perform  this  function.  Figure  6  shows  a  very  good  performance 
for  representing  a  bat  sound  signal  using  a  signal-dependent  kernel.  Radially  constant 
(rectangular)  and  radially  inverse  kernels  have  almost  as  good  resolution  and  as  a 
consequence  they  produce  some  ripples  or  artifacts  in  the  distribution.  The  radially 
Gaussian  kernel  is  a  good  trade-off  between  the  two  because  it  has  almost  the  same 
resolution  and  much  less  artifacts.  These  results  clearly  show  that  the  adaptive  optimal 
kernel  TFD  performs  much  better  than  the  fixed-kernel  methods  shown  in  Figure  5. 

The  mean-squared  error  is  used  as  quantitative  performance  evaluation  criterion.  Tables 
1  -2  summarizes  the  mean-squared  errors  between  the  TFDs  produced  by  the 
algorithms,  described  in  section  2,  and  an  ideal  distribution  (noise  free  Wigner-Ville 
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distribution).  As  expected,  the  performance  of  the  adaptive  Gaussian  chirplet 
decomposition  spectrogram  is  better  than  the  adaptive  Gaussian  expansion  spectrogram 
and  adaptive  energy  distribution.  Among  linear  TFDs,  most  of  the  window  functions 
such  as  Hanning,  Hamming,  Kaiser  and  Triangular  produce  similar  performance. 
Generally  the  bilinear  TFDs  show  better  performance  than  the  linear  transforms.  The 
tables  show  that  the  tilted  Gaussian,  the  Butterworth  and  the  generalized  exponential 
kernel  functions  produce  the  good  results  for  bilinear  TFDs.  Overall,  the  table  clearly 
shows  that  the  adaptive  kernel  approach,  especially  the  radially  adaptive  Gaussian 
kernel,  produces  the  best  results.  These  results  agree  with  our  qualitative  assessments 
(Figures  3-6). 
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Name  of  Algorithm 

Additional  Parameters 

MSE 

AED 

One  iteration 

5.4339617e-003 

AGCD 

2.2599995e-003 

AGE 

1.239271  le-002 

Linear  TFDs:  (windows  length  75) 

Rectangular  window 

2.8201 142e-002 

Triangular  window 

2.3723928e-002 

Gaussian  window 

2.8667000e-002 

Inverse  window  (14) 

p  =  0.1 

3.0140008e-002 

Inverse  window 

p  =  0.5 

3.7577926e-002 

Hanning  window 

2.4042926e-002 

Hamming  window 

2.4074127e-002 

Kaiser  window 

7  =  7T 

2.3485069e-002 

Kaiser  window 

7  =  2-7T 

2.4590848e-002 

Blackman  window 

2.6584902e-002 

Bilinear  TFDs: 

Rectangular  (constant)  kernel 

5.325361  le-003 

Bom-Jordan  kernel 

7.64486 13e-003 

Choi-Williams  kernel 

(7  =  3 

6.1759351e-003 

Choi- Williams  kernel 

CT  =  10 

4.4467070e-003 

Butterworth  kernel 

M  =  2,N  =  2,  t060  =  30 

3.4777470e-003 

Gaussian  kernel 

cr  =  100 

3.5908155e-003 

Gaussian  kernel 

a  =  200 

6.2000926e-003 

Generalized  exponential  distribution 

M  =10  ,N  =  1,  r09  o  =  30 

2.1118777e-003 

Generalized  exponential  distribution 

M  =  5,  N  =  2,  to#o  =  50 

1.6623236e-003 

Cone-shaped  kernel 

a  =  0.001 

6.8324434e-003 

Cone-shaped  kernel 

a  =  0.0001 

6.7897483e-003 

Tilted  Gaussian  distribution 

r0  =  0.1,  do  =  0.2,  r  =  0.5 

1.9228440e-003 

Tilted  Gaussian  distribution 

r0  =  0.2,  d0  =  0.1,  r  =  0.5 

8.1887901e-003 

Tilted  Butterworth  distribution 

r0  =  0.05,  60  =  0.1,  r  =  0.5 

1.131 5364e-003 

Tilted  Butterworth  distribution 

r0  =  0.1,  do  =  0.2,  r  =  0.5 

8.8955521e-004 

Table  1:  Quantitative  performance  evaluation  of  TFDs 
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Name  of  Algorithm 

Additional  Parameters 

MSE 

Adaptive  kernel:  radially  constant 

p(ip)  =  100 

1.0445602e-003 

Adaptive  kernel:  radially  constant 

<r(ip)  =  200 

6.9237386e-004 

Adaptive  kernel:  radially  Gaussian 

<r  (ip)  =  100 

8.8654181e-004 

Adaptive  kernel:  radially  Gaussian 

=  200 

6.2966003e-004 

Adaptive  kernel:  radially  inverse 

a  (ip)  =  100 

1.9123454e-003 

Adaptive  kernel:  radially  inverse 

=  200 

1.0593200e-003 

Table  2:  Quantitative  performance  evaluation  of  TFDs 
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Four  targets  -  microwave  reflectors  : 
First,  second  and  fourth  stationary,  third 
is  rotating  at  constant  angular  speed 
(DREO) 


Figure  7:  A  schematic  of  the  experiment  set-up  for  studying  an  oscillating  scatterer.  The  third  scatterer  from  the  top 
is  attached  to  a  mechanical  sine-drive.  The  other  three  scatterers  are  stationary  to  provide  reference  in  the  HRR 
profiles. 


4.2  Fusion  approach 

We  demonstrate  the  performance  of  our  new  algorithm  using  the  experimental  high 
range  resolution  radar  data.  An  experiment  was  set  up  to  compare  the  new  TFD 
approach  with  other  available  TFDs.  High  range  resolution  profiles  were  collected 
using  a  stepped  frequency  wave-form  (SFWF)  radar  mode  at  X-band  between  8.9  to  9.3 
GHz.,  i.e.,  a  synthetic  bandwidth  of  400  MHz;  the  frequency  step  size  was  1  MHz.  The 
experimental  set  up  is  shown  schematically  in  Figure  7,  which  shows  the  three 
stationary  comer  reflectors  and  an  oscillating  comer  reflector.  The  test  target  was  made 
up  of  four  comer  reflectors,  three  of  which  are  stationary  to  provide  a  geometric 
reference  and  a  contrast  to  the  distorting  shape  of  the  oscillating  reflector  in  the  HRR 
profiles. 

Figure  8  shows  the  standard  HRR  image  of  the  complex  target.  The  HRR  profile  of  an 
oscillating  scatterer  is  clearly  apparent  in  the  figure.  When  the  oscillating  scatterer  is 
stationary,  a  sharp  peak  is  observed  in  the  HRR  profile  as  expected.  This  is  shown  as 
the  top  profile  in  Figure  8.  However,  when  the  scatterer  is  in  motion,  the  spread  in  the 
spectmm  tells  us  very  little  about  the  scatterer’s  bahaviour,  as  shown  as  the  bottom 
profile  in  Figure  8. 

TFDs  immediately  provide  more  information.  Figure  9  shows  the  time-frequency 
signature  of  the  radar  signal  returned  from  the  four  scatterers,  where  the  oscillating 
curve  can  be  observed  very  well.  Figure  9  shows  the  results  of  nine  TFDs.  TFDs  used 
in  this  figure  are  Wigner-Ville  distribution,  Bom-Jordan  distribution,  Choi- Williams 
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Figure  8:  HRR  profiles  of  four-scatterer  target.  Top  profile:  all  four  scatterers  stationary;  bottom  profile:  three 
stationary  scatterers  and  one  oscillating  scatterer. 


distribution,  bilinear  TFD  with  Butterworth  kernel,  bilinear  TFD  with  2-dimensional 
Gaussian  kernel,  generalized  exponential  distribution,  cone-shaped  distribution,  tilted 
Gaussian  distribution  and  tilted  Butterworth  distribution.  These  nine  TFDs  are  used  for 
the  algorithmic  fusion  method  described  in  Section  4.  Figure  10  shows  the  resulting 
TFD  obtained  by  using  the  algorithmic  fusion  method.  An  ideal  distribution,  i.e., 
noise-free  Wigner-Ville  distribution,  is  manually  plotted  in  Figure  1 1  as  a  reference  for 
comparison.  As  expected,  Figure  10  shows  considerable  improvement  over  other  TFDs 
presented  in  Figure  9.  We  can  see  that  Figure  10  achieves  the  same  shaipness  as  the 
reference  image.  These  results  demonstrate  that  the  new  algorithmic  fusion  method 
performed  well  in  achieving  improved  resolution,  highly  concentrated  and  readable 
representation,  without  the  auto-term  distortion  and  cross-term  artifacts  that  are 
apparent  in  all  the  representations  of  Figure  9. 
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Figure  9:  Results  of  nine  TFDs:  (a)  Wigner-Ville  distribution,  (b)  Born-Jordan  distribution,  (c)  Choi-Williams 
distribution,  (d)  bilinear  TFD  with  Butterworth  kernel,  (e)  bilinear  TFD  with  2-dimensional  Gaussian  kernel,  (f) 
generalized  exponential  distribution,  (g)  cone-shaped  distribution,  tilted  Gaussian  distribution  and  (h)  tilted 
Butterworth  distribution.  The  horizontal  axis  is  time  and  the  vertical  axis  is  frequency. 
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Figure  10:  TFD  obtained  from  the  algorithmic  fusion  method.  The  horizontal  axis  is  time  and  the  vertical  axis  is 
frequency. 


Figure  11:  An  ideal  time-frequency  distribution.  The  horizontal  axis  is  time  and  the  vertical  axis  is  frequency. 
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5.  Conclusion 


This  report  evaluates  and  compares  the  qualitative  and  quantitative  performances  of 
different  time-frequency  distributions  (TFDs)  developed  for  the  past  ten  years  and 
proposes  a  new  approach  in  time-frequency  analysis  applied  to  synthetic  aperture  radar 
imaging  and  target  feature  analysis. 

Among  linear  TFDs,  most  of  the  window  functions  such  as  Flanning,  Flamming,  Kaiser 
and  Triangular  produce  similar  performance.  Generally  the  bilinear  TFDs  show  better 
performance  than  the  linear  transforms.  The  tilted  Gaussian,  the  Butterworth  and  the 
generalized  exponential  kernel  functions  produce  good  results  for  bilinear  TFDs.  The 
performance  of  the  adaptive  Gaussian  chirplet  decomposition  spectrogram  is  better 
than  the  adaptive  Gaussian  expansion  spectrogram  and  adaptive  energy  distribution. 
Overall,  the  study  clearly  shows  that  the  adaptive  kernel  approach,  especially  the 
radially  adaptive  Gaussian  kernel,  produces  the  best  results  from  all  tested  methods 
when  the  mean-squared  error  is  used  as  the  performance  evaluation  criterion.  These 
results  agree  with  our  qualitative  assessments.  A  short  review  of  different 
time-frequency  approaches  is  also  provided. 

A  new  TFD  algorithmic  fusion  method  is  presented  and  evaluated  on  the  experimental 
high  range  resolution  radar  data  and  simulated  data.  It  is  shown  that  the  TFD 
algorithmic  fusion  method  provides  an  effective  method  of  achieving  improved 
resolution,  highly  concentrated  and  readable  representation  without  the  auto-term 
distortion  and  cross-term  artifacts.  This  method  is  suitable  for  HRR  and  1SAR  data 
where  multiple  scatterers  are  present,  and  noise  along  with  artifact  reduction  are 
essential  for  target  identification  applications.  Analysis  of  the  time-varying  Doppler 
signature  in  the  joint  time-frequency  domain  can  provide  useful  information  for  target 
detection,  classification  and  recognition.  We  anticipate  that  this  new  approach  will  find 
a  wide  range  of  uses  and  will  emerge  as  a  powerful  tool  for  time-varying  spectral 
analysis. 

Our  study  shows  how  the  TFD  algorithmic  fusion  method  provides  useful  data 
visualization  of  non- stationary  signals  in  the  time-frequency  plane.  For  useful 
application  of  this  method,  however,  the  derivation  of  measures,  statistics,  or 
parameters  from  the  time-frequency  plane,  are  required  for  further  analysis  and 
processing;  for  example,  automated  target  detection  and  classification.  This  is  the 
essential  step  to  enable  TFD  algorithmic  fusion  method  to  be  employed  for  processing 
tasks  other  than  data  display.  One  of  the  advantages  of  the  TFD  algorithmic  fusion 
approach  is  that  the  noise  tends  to  spread  out  its  energy  over  the  entire  time- frequency 
domain,  while  target  signals  often  concentrate  their  energy  on  regions  with  limited  time 
intervals  and  frequency  bands.  Thus,  signals  embedded  in  noise  are  much  easier  to 
recognize  in  the  joint  time-frequency  domain.  Hence,  with  constant  false-alarm  rate 
(CFAR)  detection,  signals  can  be  detected  and  reconstructed  by  using  only  the  detected 
time-frequency  coefficients.  Therefore,  the  target  signal  buried  in  noise  can  be  detected 
and  its  parameters  can  be  measured.  By  applying  time- varying  frequency  filtering,  the 
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signal-to-noise  ratio  (SNR)  can  also  be  enhanced  and  the  resulting  statistics  would 
outperform  the  conventional  measure.  We  intend  to  investigate  these  scenarios  in  the 
future  study. 

Another  significant  advantage  of  the  TFD  algorithmic  fusion  method,  as  we  have 
already  seen,  is  that  it  allows  us  to  determine  whether,  or  not,  a  signal  is 
multicomponent  or  not  and  its  ability  to  decompose  the  signal  in  the  time- frequency 
plane.  That  is,  the  time-frequency  analysis  enables  us  to  classify  signals  with  a 
considerably  greater  reflection  of  the  physical  situation  than  can  be  achieved  by  the 
Fourier  spectrum  alone.  Based  on  this  time- frequency  localization,  one  may  wish  to 
select  “desired”  signal  components  and  remove  unwanted  contributions,  by  applying  an 
appropriate  mask  function.  Time-varying  filtering  and  signal  synthesis  have  been 
applied  to  many  signal  processing  areas  and  we  intend  to  investigate  these  approaches 
in  our  future  study. 

This  work  is  especially  relevant  to  the  FIRR/ISAR  imaging  capability  of  the 
surveillance  radar  systems  on-board  of  the  CF  C-140  Aurora  patrol  aircraft  and  the  US 
Navy’s  P3  Orion  surveillance  aircraft.  Moreover,  target  recognition  based  on  radar 
imagery  will  play  an  active  role  in  future  CF  initiatives  on  1SR  (Intelligence, 
Surveillance  and  Reconnaissance)  for  land,  air  and  maritime  applications. 
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